Computing the Bézier Control Points of the Lagrangian Interpolant in Arbitrary Dimension

The Bernstein-Bézier form of a polynomial is widely used in the fields of computer aided geometric design, spline approximation theory and, more recently, for high order finite element methods for the solution of partial differential equations. However, if one wishes to compute the classical Lagrange interpolant relative to the Bernstein basis, then the resulting Bernstein-Vandermonde matrix is found to be highly ill-conditioned. In the univariate case of degree n, Marco and Martinez [18] showed that using Neville elimination to solve the system exploits the total positivity of the Bernstein basis and results in an O(n2) complexity algorithm. Remarkable as it may be, the Marco-Martinez algorithm has some drawbacks: The derivation of the algorithm is quite technical; the interplay between the ideas of total positivity and Neville elimination are not part of the standard armoury of many non-specialists; and, the algorithm is strongly associated to the univariate case. The present work addresses these issues. An alternative algorithm for the inversion of the univariate Bernstein-Vandermonde system is presented that has: The same complexity as the MarcoMartinez algorithm and whose stability does not seem to be in any way inferior; a simple derivation using only the basic theory of Lagrange interpolation (at least in the univariate case); and, a natural generalisation to the multivariate case.